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^ _ ISOTROPY CONDITIONS FOR LATTICE BOLTZMANN SCHEMES. 

■ APPLICATION TO D2Q9 
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Abstract. In this paper, we recall the linear version of the lattice Boltzmann schemes in the frame- 

^-rfH ■ work proposed by d'Humieres. According to the equivalent equations we introduce a definition for a 

1 J* , scheme to be isotropic at some order. This definition is chosen such that the equivalent equations are 
^— ( ' preserved by orthogonal transformations of the frame. The property of isotropy can be read through 
(-H I a group operation and then implies a sequence of relations on relaxation times and equilibrium states 

. that characterizes a lattice Boltzmann scheme. We propose a method to select the parameters of the 

^ ' scheme according to the desired order of isotropy. Applying it to the D2Q9 scheme yields the classical 
constraints for the first and second orders and some non classical for the third and fourth orders. 
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Introduction 



Numerous phenomena, in particular acoustic properties, are modelled by partial differential equations (PDEs) 
like Navier-Stokes equations. One of the main property of these PDEs is to take into account the isotropy of 
the environment. Namely, if the physical properties of the environment do not depend on the orientation, 
I ' the associated PDEs are invariant by orthogonal transformations of the spatial frame. In this paper, we are 
I interested by lattice Boltzmann schemes [5l[6l[9] used by the industrialists because of their efficiency. As the 
1^ ' directions of the lattice are privileged, they are not a priori numerically isotropic. Because of the criss-cross 
. ^ ] pattern of the scheme DdQq (where d represents the space dimension and q the number of discrete velocities), 
■ it is natural to look for parameters that give the same behaviour on each axes (i.e by a rotation of angle 7r/2 in 
^ , 2D for example) and by extent we want to obtain this result for each orthogonal transformation of the spatial 
- - ' frame. Our purpose is then to understand and control the lack of isotropy for these schemes. 

The first thing to do consists in defining the isotropy for lattice Boltzmann schemes. Actually, there is no 
uniform way to define isotropy: for example in [3] the authors focus on the maintain of isotropy of the stress 
tensor required by fluid equations. An other point of view is described in [3]: an isotropic fluid is simulated by 
eliminating as much as possible the effects due to the symmetry of the lattice and it means solving an eigenvalue 
problem. The problem of the space transformation is evoked in [TOj[Tl] : in [11] rotating frames are considered 
and in |10| the authors introduce a deflnition of isotropic diffusion that constraints the collision operator to be 
invariant under all isometrics of the Bravais lattice. In this paper, we deflne isotropy as the invariance by the 
orthogonal transformations of the frame and we express it on the equivalent equations. Namely, the equivalent 
equations are PDEs that come from a Taylor expansion of the lattice Boltzmann scheme [5] and a PDE is said 
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isotropic if it is invariant by rewriting it in the "new" frame. Finally, a DdQq frame is isotropic at order M if 
its set of equivalent equations of order n, 1 < n < M , is isotropic. 

In the first Section we quickly recall the lattice Boltzmann method and give some results on the equivalent 
equations. 

In the second Section, we exploit the definition of isotropy for PDEs in an algebraic point of view. More 
precisely, we first define a new frame by orthogonal transformation. After redefining all the objects in the new 
frame, we identify both of the equivalent equations: the old one and the new one. The isotropy is guaranteed 
as soon as the various tensors in the PDEs appear as a fixed point of some group operation. Then, we obtain 
a very useful set of relations on the coefiicients of the equivalent equations. 

In the last Section, a lattice Boltzmann scheme is said isotropic at order M if its set of equivalent equations 
of order n, 1 < n < M, is isotropic in the sense explained in Section [5] Since the coefiicients of the equivalent 
equations depend on the parameters of the scheme (relaxation times and equilibrium states), we are able to 
characterize, for the related scheme, the lack of isotropy on these parameters for each order. More precisely, 
we are able to determine order after order the relations between the parameters of the scheme that have to be 
satisfied to have isotropic equivalent equations. We determine these relations for the D2Q9 scheme and specify 
all the sets of parameters that give isotropy until the fourth order. 

The lattice Boltzmann schemes considered in this paper are linear (for applications reasons) while the pro- 
posed definition of isotropy is more general and can be used as soon as the equivalent equations of the scheme 
are known. 



1. Lattice Boltzmann method and equivalent equations 

In this section, we recall some notations related to the lattice Boltzmann method. Then, we introduce the 
equivalent equations as well as their link with the physic. This link is essential because it allows us to decide 
what parameters have physical values and what parameters have to be chosen in a allowed set of parameters. 
We will see below (in section[3]) what for a choice we have for this set of parameters in order to guaranty isotropy 
of the scheme. 



1.1. Notations and description of the lattice Boltzmann method 

Wc use the scheme proposed by d'Humieres in [5]: we consider a regular lattice C with typical mesh size 
Ax. The time step A< is determined thanks to the velocity scale A: At ~ Ax/X. For the DdQq scheme, we 
note V = {vj)i<j<q the set of q velocities and we assume that for each node x of C, and each vj in V, the point 
X + VjAt is also a node of the lattice C. The aim of the DdQq scheme is to precise the particle distribution 
/ = ifji^T^))\<j<q (where * represents the transpose) for a; G £ and discrete values of time t by solving a 
discretization of the Boltzmann's equation in two steps for every time: collision and transport. We finally have 
to propagate the so-called lattice Boltzmann scheme: fj{x + VjAt,t + At) = f*{x,t), V.t G C, 1 < j < q, 
where /* represents the jth component of the particle's distribution after the collision. In order to implement 
this scheme, as it is done in [5], we introduce the orthogonalized moments m defined by m := M f where M is 
a given and invertible matrix. Concerning the moments m, we distinguish those who are preserved during the 
collision, called W 6 M^, of the complementary ones denoted by Y" S M''"^. Since we are interested in linear 
applications, we assume that the collision operator is linear, so there exists a matrix J such that m* := Jm, 
with the following structure: 

( W* \ , /Wat \ f W 

m = { ^ = Jm = ^^j^ Id,_^-5 ) \ Y 

where E G Aiq-N.N{^) and 5* is the diagonal matrix of the relaxation times Sk = At/vk, N + 1 < k < q, with 
< Sfc < 2, + 1 < < g for stability conditions. We introduce a last classical notation: = 1/sfc — 1/2 e 
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R+, 1 < k < q. Finally a time step of linear lattice Boltzmann scheme scheme reduces to: 

f{x + v^At,t + At) ^ Ar^ J M f{x,t), VxeC, I < j < q. (1) 

Thus, the parameters of the lattice Boltzmann schemes are the coefficients of the matrices S and E. This 
coefhcients are usually chosen thanks to the experiment (by the physicists and the industrialists). The aim of 
this paper is to precise these coefficients according to the property of isotropy. 

Remark 1.1. There exist three types of lattice Boltzmann schemes according to the values of the relaxation 
times: either all of relaxation times are equal (BGK for Bhatnagar - Gross - Krook scheme) or there is only two 
different relaxation times (TRT for Two Relaxation Times) or there is no conditions on these (MKT scheme for 
Multiple Relaxation Times). Since we investigate the behavior of the scheme in the sense of the isotropy, it is 
natural to consider the third case: it is the most general case even if it is not the most often used. Thus, there 
are q — N independent unknowns parameters on S. 

Focus now on the matrix E that is chosen full in order to give us as degree of freedom as possible: it gives 
N{q — N) additional parameters. Actually, it is not an usual choice: the matrix E is usually chosen sparse by 
physicists, but in our knowledge, there is no mathematical justification to this type of assumption (based on 
the physical "good sense"). It is remarkable, that the isotropy condition for the DdQg scheme can give some 
justification of this set of coefficients. 

Finally, both matrices E and S give {N + 1) x [q — N) parameters to be precise in order to conserve physicals 
and isotropic properties. 

For example, for the D2Q9 scheme, we consider the classical moments: the density p and the momentum 
Q i^xTly) that are conserved during the collision and the energy e, the square of the energy £2, the flux 
{(fxT^Py)^, and the tensor {ipxx,'^xyT ■ Then, we orthogonalize them like it is done in [5]. For this choice of 
moments we define both equilibrium states and relaxation times matrices by: 

/ EPX^ EP.X EP^X^ EP^X^ EP^^X' EP^^X^ \ 
E^ := EI^X EI^^X^ E^iX^ E^j^^X^ E^i^X E^i^X and 5 
V E^'X EtlX"^ E^X^ E$IX^ E$l^X E^^X ) 

1.2. Equivalent equations 

We first introduce the notation for the element of • that is on the intersection of the ith line and the 
jth column. By convention, a Latin letter is an index related on the moments (of size N) while Greek letter 
represents an index related on the space dimension (of size d). We then introduce the system of N equivalent 
equations [6l[7] of order M with which the linear lattice Boltzmann scheme is consistent at order M : 

M 

dtW + ^nA' : V^W ^ 0, l<i<N, (2) 

n=l 

where ^ < n < A/, 1 < i < are tensors of order n that take into account the coefficients of E and S. 

By convention, the maximal contraction operator : is defined by: 

(„A : V"W)* := : V"Ty := ^ „4"^'- '""aa^ • • • 9a„ W^ l<i<N,l<n<M. 

l<j<N 
l<ai ,■ ■ ■ ^ctji <d 

These equivalent equations come from a formal calculus explained in [8] and an algorithm that is easy to use 
is described in [5]. 

Remark 1.2. The tensors „A belongs to the space M^^"'', 1 < n < M, but because of the Schwarz property, 
numerous of coefficients of are equal (or are proportional between them). In order to count how much 



■ -Di^SC-^ei "^£2' "^tpxl ^<fiy^ "^tpxa:' ^ipxy)' 
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of these coefficients are independent, we have to calculate the number of d-uplet • • • , (3d) G N'' such that 
= n. This is a well-known combinatorics result: this number is the number of d combinations with 
repetition of a set which cardinal number is n and it reads: FJJ := {d + n — 1)!/ {nl{d — 1)!) and „ A G iS'(7i, d), 
where S{n, d) is the space of the symmetrical tensors of " . 

The equivalent equations ([2]) are also very useful to make the link with the physics. For example, it is well- 
known (see for example |7|) that there exists a D2Q9 scheme {i.e. a choice of moments S M^) such that Q 
is consistent with Navicr-Stokes equations at the order 2. Then, there exists a relation between both viscosities 
(shear viscosity ( and bulk viscosity /i) and two relaxation times (see for example |5]). Furthermore, equilibrium 
state of the energy given by EPp+ E^^Qx + E^^ qy + Ej^q^ is such that E^ only depend on the sound velocity. 

Finally, our purpose is to precise as much as possible the (iV -I- l)(g — A^) — 3 parameters of the lattice 
Boltzmann scheme in order to understand the lack of isotropy of these schemes: for example, in the D2Q9 
scheme, parameters that are given thanks to the environment are E^ (it depends on the sound velocity), (it 
depends on the bulk viscosity C,) and ct,^^^ (it depends on the shear viscosity ji) that means 21 parameters. 

2. Algebra for rotating the equivalent equations 

We first recall the definition of an isotropic PDE: 
Definition 2.1. A PDE 

M 

dtW + nA' : V'W = 0, 1 < i < A^, (3) 

n=l 

is said isotropic if it is invariant by orthogonal transformation of the frame. 

We introduce the matrix r G C'd(R) that represents an orthogonal transformation of the frame (O, e). We 
then use the notation • for every clement in the "new" frame: in this way, if the vectors a; G R'^ and the 
moments W G are known in the frame (O, e), we have to define the news moments W, the new vector 
X (such that x := rx), the new basis e, the news tensors „A*, 1 < 7i < M, I < i < N and the new 
partial derivatives da, 1 < a < d. Finally, let be R{r) the matrix in C'jv(R) such that W :~ R{r)W. 
More precisely, R : C'd(R) C?jv(K) defines a morphism that gives a group representation. In the present 
contribution, mass and momentum are conserved thus W = {p, q)^ , N = d+1 and for every r G ©^(M), we get 
R{r) = Diag(l,r) G OAr(R). 

The aim of this section consists in giving the relation between objects that are defined before and after the 
transformation of the frame. First, we give a relation between both partial derivatives and we deduce a relation 
on both derivatives of the moments. Second, we introduce a group operation via an homomorphism that will 
show itself very useful for the definition of the isotropy. 

Lemma 2.2. For every 1 < a < d and I < i < N , the derivatives daW at the ith component of the moment 
W reads: 

- daiRir)Wy = ^i?}(r)(r-i)^^Ty'. (4) 

Proof. As a; = r 5, for each 1 < a < d, we have the following relation on the partial derivatives da along the 
ath direction: 

l</3<<i i<P<d 

Then, it is easy to characterize the partial derivatives of the moments and we obtain the desired result (|4]). □ 



We also define a group operation through the following Proposition. 
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Definition - Proposition 2.3. Let n be in N* and the map $„ : Od{^) x S{n,d) S{n,d) defined for 
1 ^ *) i N , 1 < ai, • • • ,an<d, 1 < n < M , by the relation: 

($„(r)(„A));.'—"" := (5) 

l</3i,--- ,/3„<d 
l<k,l<N 

It is an homomorphism that characterizes a group operation. 

Proof. Let Id be the identity matrix of 0^(1^), tlien for every „A e S{n, d) we have <I>„(Id)(„y4) = „A. It remains 
to prove the following property: Vr, r' e C'd(R), V„A £ S{n,d), $„(rr')(„A) = $„(r) ($„(r')(„A)). First, by 
using the definition of the homomorphism $„ given in ([5]), we get for 1 <i,j < N and 1 < ai, • • • , a„ < d: 

l</?i,---,ft.,7ir--,7„<<i 
l<fc,i,m,o<Ar 

Since i?(rr-') = R{r)R{r'), and so i?(rr')^^ = R{r')-'^R{r)-\ ^ gives: $„(r) ($„(r')(„A)) = $„(rr')(„A) 
and the homomorphism (/)„ truly defines a group operation. □ 

Property 2.4. Let ?' be an orthogonal transformation of Od{^), then the system ([3]) reads in the new frame 

dtW' + nA' : V'W = 0, 1 < i < iV, (7) 

n=l 

where „A = <I>„(r)(„A), 1 < n < M, or equivalently „A = <I>„(r~^)(„A), 1 < n < M. 

Proof. Let us write ([3]) after the orthogonal change of frame r. Because of the equality W = R{r)W , dtW 
becomes dtR{r) W. Using the property (|H), we obtain a relation that gives : V"W^ in function of the partial 
derivatives of the moments in the new frame, more precisely for 1 < i < iV we get: 

l</3i.---,/3„,7i,---,7„<d 
l<k,l.m,o<N 

Thus, by identifying both ^ and R{r) x (O, we characterize „A, 1 < Ji < M , thanks 

„I*:V"W-= Y : V"Ty) , VVT eR^, 1 < n < A/, 1 < I < iV, 

and that gives the desired result. □ 

Finally, we obtain a characterization for a PDE to be isotropic: 

Corollary 2.5. Let ^ be a PDE of order M, it is isotropic if and only if nA = nA, l<n<ALVre 0d(R) 
and it means that the tensor nA is a fixed point of $„, 1 < < that is $„(r)(„A) = „A, Vr e Od{^), 
l<n< M. 
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3. ISOTROPY CONDITIONS FOR LATTICE BOLTZMANN SCHEMES 

In this section, we first give a definition of isotropy for lattice Boltzmann schemes and then we give an 
application on the D2Q9 scheme. 

Definition 3.1. A lattice Boltzmann scheme is said isotropic at the order M if the system of equivalent 
equations ^ at the order M is isotropic. 

Furthermore, we denote by L^lr) := X]i<n<Af (*^n(^)(n^) ^ the lack of isotropy at order N of the 

orthogonal transformation r of the equivalent equations ([2|). 

Then, we use this definition on the D2Q9 scheme in order to obtain all the sets of parameters that improve 
this scheme (in the sense of the isotropy) until the fourth order. In this case, the orthogonal transformation r 
is the well-known rotation matrix, d is equal to 2, iV = 3 and M = 4. 

Proposition 3.2. Let L^ir) be the lack of isotropy for the D2Q9 scheme at fourth order for the orthogonal 
transformation r, then we get: 

• U{r) = 0(At2),Vr G iff E^^=^ = E^^ = E^,^^ = i??^^, - - E'^^^ = E^i^ = E^l^ = 0. 
. U{r) = 0(At3),Vr e O^i^) tff EP^ = E^ = = i??^^ = and E^i = E^^. 

• Ulr) = 0(At4), Vr e 02(R) iff (E^^^ = E^ = 0, = a^^^, E^i = -I) and (either 2EP^ +A + 3EP = 
ora^^ = a^^ ^ l/(12cr^^Jj. 

• L^r) = 0(At5), Vr G e>2(M) iff (2EP^ + 4 + 3EP = 0, a, = a^^ = a^^ = ^/{6a^^J) and (either 
2 + SEP ~ or (Te2 = cfip^-^). Thanks to this last condition gives the following relation on both 
viscosities: C = 6/i. 

Remark 3.3. Some of the properties for an isotropic D2Q9 scheme that are obtained here are well-known and 
naturally used in this scheme. In fact, the parameters that are taken null in order to obtain isotropy at first 
and second orders are also taken null for isotropic reasons considering the kinetic solution at equilibrium in the 
continuous environment. However, the results on the third and fourth orders are more surprising even if some 
usual D2Q9 schemes can be seen as particular cases of those we propose. 



The proof of Proposition 13.21 is obtained thanks to a rigorous implementation of isotropy conditions with a 
software of formal calculus and [3] give details of it. 
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